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The quasistatic limit of the antiplane shear-wave speed ('effective speed') c in 2D periodic lattices 
is studied. Two new closed-form estimates of c are derived by employing two different analytical 
approaches. The first proceeds from a standard background of the plane wave expansion (PWE). 
The second is a new approach, which resides in x-space and centers on the monodromy matrix 
(MM) introduced in the 2D case as the multiplicative integral, taken in one coordinate, of a matrix 
with components being the operators with respect to the other coordinate. On the numerical side, 
an efficient PWE-based scheme for computing c is proposed and implemented. The analytical and 
numerical findings are applied to several examples of 2D square lattices with two and three high- 
contrast components, for which the new PWE and MM estimates are compared with the numerical 
data and with some known approximations. It is demonstrated that the PWE estimate is most 
efficient in the case of densely packed stiff inclusions, especially when they form a symmetric lattice, 
while in general it is the MM estimate that provides the best overall fitting accuracy. 



I. INTRODUCTION 

Effective material properties of composites have been 
and remain a topic of much interest in micromechanics, 
see the reviews!®]. The recent surge of research into 
the properties of metamaterials and phononic crystals 
has heightened attention, particularly for periodic sys- 
tems. In this context, considerable work has been done 
on the low-frequency, or quasistatic, limit of the antiplane 
shear-wave speed in 2D periodic structures (referred to 
as the 'effective speed' c in the following; note that this 
value also yields the limit of the fundamental velocity 
branch of shear plate waves). A natural tool for tackling 
the problems with periodicity is the plane-wave expan- 
sion (PWE) . An explicit PWE-expression of the effective 
speed c via an infinite sum of Fourier coefficients has been 
obtained ire' and was broadly used afterwards for com- 
puting c in various periodic materials. Note some other 
semi-analytical techniques that were used for numerical 
evaluation of c, such as scaling^ and mixed-variationaP 
methods. In turn, the multiple-scattering theory (MST), 
which deals directly with the inclusion/matrix bound- 
ary problem, has proved to be expedient for deriving the 
effective speed in an approximate but closed form. By 
means of MST, such simple (in appearance) estimate of 
c, which had been known for certain statistically uni- 
form models in micromechanics, was recently extended 
to phononic crystals with a periodic microstructure of 
inclusionpHini 

The main results of the present paper are concerned 
with both analytical and numerical aspects of the prob- 
lem of evaluating c. The analytical development aims at 
finding approximations of c by means independent of the 
MST. The starting point is the general expression for c in 
the operator form that may be further specialized to ei- 
ther Fourier space or x-space. On this basis, we provide 
two new closed-form estimates of c derived by employ- 
ing two different analytical frameworks. The first is the 
PWE approach, which commences from the formula of 4 . 



The second is a completely new approach based on the 
monodromy matrix (MM), which is a fundamental ob- 
ject for the ID periodic problems (cf. the state-vector 
formalism) but it has not seen much, if any, far-reaching 
application in 2D. Here the MM is introduced as a multi- 
plicative integral, taken with respect to one coordinate, of 
the matrix with components defined as operators acting 
on the functions of the other coordinate. On the numer- 
ical side, we develop an efficient PWE-based scheme for 
computing c, in which the matrix inversion is replaced by 
the power series that is judiciously gauged for its faster 
convergence. The results are applied to several examples 
of 2D square lattices consisting of two and three high- 
contrast components with filling fractions /, for which 
the new PWE and MM estimates of c (/) and its known 
MST estimate are compared against the benchmark of 
the numerically computed c (/) . In brief, it is demon- 
strated that the PWE estimate is efficient in the case of 
densely packed stiff inclusions (where the MST estimate 
fails) and is particularly useful for the symmetric binary 
lattices invariant to interchanging their components (in 
which case the MST formula is ambiguous); but it is the 
MM estimate that provides the best overall fit over vari- 
ous lattice configurations considered. 

The paper is organized as follows. The background 
expression for c is presented in EjTTJ The PWE and MM 
closed-form estimates of c are derived in £ III The numer- 



ical scheme used for computing c is described in §IV| Ap- 
plication of analytical and numerical results to 2D square 
lattices is discussed in fV) Concluding remarks are pre- 
sented in §VI| Appendix expands on the convergence of 
the implemented numerical scheme. 



II. BACKGROUND: EXACT EXPRESSION FOR 
EFFECTIVE SPEED 

Consider a 2D periodic locally isotropic medium with 
the density p (x) and the shear coefficient ^t(x), which are 



2 



real positive piecewise continuous functions satisfying 

,2 



p[x- 
p(x- 



n j a j) = K X ), 



/i(x) 



(1) 



for any x £ R 2 , rij € Z and some linear independent 
translation vectors a.j € R 2 that form the irreducible unit 

cell T = Ylj=i tj a j (.tj £ IP) 1]) °f the 2D periodic lattice. 
Let ei,e2 = {e^} be an orthonormal base in R 2 , and 

x • y = Yli=i x iVi be the scalar product in R 2 , where Xi 
are the coordinates of an arbitrary vector x with respect 
to {ej}. Denote 



where a 7 • 



= (A- 1 ) i e,, g = 
6 jk (j,k = 1,2), 



X) i=1 27m i b i. ( 2 ) 



means transpose, 



and g = Ylj—i gj&j is the reciprocal lattice vector whose 
components (51,32) in { e j} take all values from the set 
r = 2tt (A _1 ) T Z 2 . In the following, the Fourier coeffi- 
cients of a periodic function /(x) are indicated by a hat: 



/(x) = J2 /(g)e iS ' X 



/(g) 



1 

ITI 



/(x)e-* x dx=(/(x)e-^ x >; (3) 



and the same notation (•, •) is used for the scalar products 
in g-space and in L 2 (T): 



(f,h) 
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/(g)M-g) 



^Jj(x)h* (x)dx 



(A*)- (4) 



The antiplane time-harmonic displacement t; (x, t) = 
v (x) e _iLJt is determined by the wave equation 

V • (^i(x)Vw(x)) = -p(x)uj 2 v(x) (5) 

with periodic p(x) and p(x). By this periodicity, v(x) = 
w(x)e lk ' x where u (x) is periodic and k = fere (|re| = 1) is 
the Floquet vector, and so Eq. ([5| can be cast as 

(Co + C\ + C2)u = pui 2 u with Cqu = — V(/iVu), 
C\u = — ik • (/iVit + V(/xu)), C2U = k 2 pu. (6) 

To find the effective speed c(re) = lim w fc-^o w (k) /fc, con- 
sider the asymptotics u 2 = uj 2 + ui 2 + u 2 + 0(k 3 ). It 
is evident that w 2 , = is an eigenvalue of (JgJ) with the 
eigenvector uq = 1. Therefore perturbation theory yields 



U17 



(Ciito,u ) 
(pu ,u ) ' 



(C 2 w ,w ) - (C 1 C 1 u ,C 1 u ) 
(pu ,u ) 



(7) 

where (puo,u Q ) = (p) and (£2^0, w ) = k 2 (p). Note that 
(CiUq, uq) — by periodicity of /i, hence w 2 = and the 



operator Cq 1 is defined in the subspace L\ of all functions 
/ orthogonal to 1, i.e., such that (/) = 0. Thus c 2 (re) is 
expressed via the averaged density (p) (= Jo (0)) and the 
effective shear coefficient /u e ff(re) as follows: 

c 2 (k) = McffW /ieff(K) = (M) _ M(K) 

^ g)*-* <»> 

The operator Cq 1 is compact and also self-adjoint and 
positive, whence 



c 2 ( K )<(p)/(p) 



(9) 



It is worth emphasizing that the perturbation theory 
enables an efficient shortcut to an explicit expression of 
the effective speed c, in which the quadratic form M (re) 
may be specialized to either g- or x-space. Taking a 

double Fourier expansion of mo, 



M(k) = ]T M(g)/i(-gO(g-K)(g'-K) 
^ — 'g,g'er\{o} 



x (p (g - g') S • g') 1 



(10) 



provides the PWE-representation of c 2 (re) obtained irP. 
Viewing Eq. ([8| along with the equation Coh — dp/dxi in 
x-space is precisely equivalent to the formulation of qua- 
sistatic limit by the scaling approach, see^R The above 
derivation, taking a few lines, does not need the scaling 
ansatz. Moreover, while the central point of the scaling 
approach is the use of the Fredholm alternative (Lemma 
1 in Ch. 4 ofS), the same is inherent to the perturba- 
tion theory 'by construction' whereby the eigenfunction 
perturbations are confined via the operator Cq 1 to the 
subspace L\ orthogonal to the unperturbed eigenfunc- 
tion (in g-space, this is implied by the summation over 
g € T\{0} in Eq. flit)])). Finally, note that Eq. fls) can 
be further developed by using the monodromy-matrix ap- 
proach, see §111 B| 



III. ESTIMATES OF THE EFFECTIVE SPEED 



PWE estimate 



Eq. ( 10 ) of 4 defines M (re) as a scalar product in the 
Fourier space I 2 (T\ {0}), 



M(re) 



(B x d,d) 



(11) 



where B is an infinite matrix and d an infinite vector 
with components 

B=(B [g, g']) g , g , er \{0} : B [S, g'] = ? (g - g') g ' g'; 

d=(d(g)) : d(g)=/2(g)g-re. (12) 



By definition (12[), B 1 is a compact operator in 
I 2 (r\ {0}) . Let us further cast /i(x) in the form 



p(x) = p + A*a(x) 



(13) 
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where /io is some positive constant and hence 
M (g - g') = Mo<W + MA (g - g')- Denote 

C(^ ) = (C [g,g']) g!g , erU0} : 

C[g,g 

D = diag(|g|) gerx{0} ; 
f = D _1 d = f /(g)) 



MA/ /^_g_ _g_ 
Mo lg gJ |g|'|g'| 



/(g) = M (g) |fr • « = Ma (g) jfr • k; 

|g| |g| 

= Y~ 2 F v k,k, (F a =F„) (14) 
It follows from ([II]), ([l2j> and ([13]), Q that 
B = Mo D (I + C) D, M(k) = Mo" 1 ((I + C)- 1 f , f) , 

where I is an infinite identity matrix. Note that I + C is 
positive and that it satisfies the identities 



(i + cr 1 = v m (-c) n + (-c) m+i (i + cr 1 

* — ^n— 

+ ((-Cr+^I + CJ^f.f 



(16) 



Taking (16>) with m — yields 



M(«) = Mo ^(«) - Mo 1 (C (I + C)- 1 f , f 



(17) 



Consider ( 171 for two different choices of /io > 0. If /Jo = 



max /<(x) = /i max then /iA (x) = m( x ) — Mo is negative, 
hence so is C and therefore the second term on the r.h.s. 
of |l7| ) is positive. If /io = min ji (x) = Mmin then the 
above signs are inverted. Thus cL^( K ) — M(k) < 
Pmln^( K )- Combining this with dst ) gives the bounds 

(m) - < Mcff(K) < (m> - for any /i(x). (18) 

Mmin Mmax 

The lower bound is not very interesting since it may be- 
come negative if /i m i n is small. The upper bound rein- 
forces the inequality ^ as 



2/ s 1 ( , , -F(«0 
C 2 M < y-r M - 

(M) V Mmax 



(19) 



It is natural to inquire as to what choice of /io provides 
the best estimate of Meff( K ) within the bounds ( fl8| ). To 
answer this question, let us formally consider Eqs. ( 16 ) 
truncated as follows: 



(i + cr 1 ^™ (-C) B , 
MW^^" 1 ((-C)"f,f), 

*—^n=0 



(20) 



The sufficient condition for convergence of both series as 
77i — ^ oo is ||C|| < 1, where ||-|| is an operator norm. 
Hence we need to take /i which minimizes ||C(/io)||. 
Note from ( 14 ) that C is close to the operator of mul- 
tiplication by /i a (x) //io , so ||C|| may be gauged by the 
value max x |/iA (x) /mo|- Its minimum over all choices of 
/i is reached when Mo — I (Mmax + Mmin)- Thus a simple 
estimate, given by a single first term M(k) ~ F(k)/(1o of 
(20>), can be taken as 



c 2 (k) 



1 



M (ftp 

(p) ~ (P> 



(M) 



Mo 



with /i 



Mn 



(21) 



Note that the obtained estimation is a general result in 
the sense of having the same form for an arbitrary pe- 
riodic dependence /i (x) , but it certainly provides a dif- 
ferent accuracy for different types of m( x )- For instance, 
consider two extreme examples: a stiff composite with 
small admixture of a highly contrasting soft ingredient 
and the inverse case where these two components form 
a soft material with a stiff reinforcement. The common 
ratio of geometrical progression (20>) with /io = ~p has a 



similar absolute value (gauged by max x |/iA (x) /m|) for 
both cases but is likely to differ in sign, since C is close to 
multiplying by /iA (x) = /i (x) — Jl and hence should be 
positive (negative) definite when the stiff (respectively, 
soft) component is volume dominant. Obviously a sign- 
alternating progression converges faster. Thus the PWE 
estimate, which is the leading-order term of (20 2), is ex- 
pected to be more accurate in the former case of a pre- 
dominantly stiff composite with a small volume fraction 
of a soft material and less accurate in the latter, inverse, 
case. This observation is illuminated by the examples in 



£ V A2 



It remains to supply the closed-form relations for F (k) 
From its definition in (14 1, 



trace (Fy) = £ 



g er\{o} 



lM(g) 



= ((m-(m)) 2 ) = (m 2 Mm) 2 . (22) 



Hence by (19) and (21 1 the sum of squared effective 



speeds along any pair of unit orthogonal vectors Hi in 
M 2 satisfies 



i=i 

2 



i=i 



TT 2 </i) - 2 

(P) \ Mmax 



Mmax 

<M 2 ) - (pY 



Mn 



(23) 



The quadratic form c 2 (k) is known to be independent 
of the orientation of K in M 2 if m( x ) (thus also /i (g) 
and c 2 (k)) is invariant under three- or fourfold rota- 
tions about the axis normal to the x-plane. In this case, 



4 



c 2 (k) — \ Y^i=i ° 2 ( K i) f° r an y K an d thus (23) gives 



notation 



(p) 



(p) 



w- 



(p 2 ) (p) 2 \ 



Pm&x ~t~ Pmin 



= c pwei (24) 



where the notation Cp WE is introduced for future use to 
distinguish this estimate from those obtained by other 
methods. For a piecewise homogeneous periodic mate- 
rial consisting of J = 1,2, ... components with constant 
fij, pj and with filling fractions fj (J2 fj = 1), E Q- 
(24) obviously specializes by setting (•) —J2j(')jfj an d 
Mmax/min = (max/ mm) 3 p,j. 

The above results are formulated for the 2D periodic 
media; however, they can be readily adapted for equa- 
tions similar to ([5} with x e R d of any dimension d > 2, 
e.g., for 3D equations of heat conduction or fl uid acous- 
tics. Indeed, replacing Ym=i by Yli=i keeps (|22| intact 
and replaces the factor 2 by d before (p) in ( |23p , which 
leads to c 2 (k) = \ Y^i=i ° 2 ( K i) ^ ° 2 ( K ) ^ s independent 
of k € W 1 . This is the case for d — 3 under cubic 
symmetry!^ For example, consider a 3D-periodic fluid- 
like cubic structure with bulk modulus K(x) and density 
p(x). Based on the standard equivalence between SH — > 
acoustics under the interchange p — > ^ftT 1 , p — > p^ 1 , the 
PWE bound and estimate of the effective acoustic speed 
follow in the form 



3 /3max + P m f n 



B. MM approach and the estimate 

In this subsection we develop the x-space approach 
basing on the monodromy matrix (MM). The idea implies 
casting the wave equation in matrix form containing an 
ordinary differential operator with quasi-periodic bound- 
ary condition in one coordinate, integrating this system 
using the multiplicative integral in the other coordinate, 
and applying perturbation theory to express the result 
via the scalar product in I? (T) that enable eliminating 
the operators and yields the closed-form approximate so- 
lution in the form of double integrals of /i(x) and p(x). 
Thus the MM approach is performed in x-space. 

It is convenient to assume for the moment that the 
functions p(x) and p(x) in the wave equation ([5} are 
smooth functions, which are periodic on the 2D rectan- 
gular lattice with the unit cell T 3 x = (2:1, £2) formed 
by the translations ai 2 || ei 2 (see pTl). Alongside the 



— |T| JT 



jry J T ■ dx introduced in |3} , denote 

H<-U, 2 = 0) (26) 



• da;,- 



and let, for brevity, ai^ be of unit length so that 
T = [0, l] 2 . Imposing the Floquet quasi-periodic con- 
dition along one of the coordinates, say xi, leads to 
v (x) = w (•, 2:2) e lklXl where w (-,£2) = w (xi) for any 
fixed x 2 and w (x\ ) is an absolutely continuous periodic 
function: 

w(xi) e W = {w(xi) € AC[0, 1] : 

w(0)=w (1) , w' (0) = w' (1)} (27) 

with ' meaning d/dxi. On these grounds, Eq. ^ can be 
rewritten in the form 



d 

On = n — Vi Q 
0x2 



M - X (X) 

A-uj 2 P {x) 

w(-,x 2 ) 
p(x)dw (•, x 2 ) /dx 2 



(28) 



where the operator A = A (fci, x 2 ) acting on the compo- 
nents of rj as on functions of x\ is defined in the space W 
by the definition 

A(k 1 ,x 2 )w(xx) 

= -e~ ik ^ (p{x 1 ,-){e ik ^w{x 1 ))')' 

= — ipw')' — iki (ypw 1 + (pytv)') + kfpw. (29) 

The solution rj (x) of Eq. ( [28] ) with the initial condition 
?7(xi,0) = 770(^1) at x 2 — can be represented in the 
form 

77 (x) = M [x 2 , 0] 770 (xi) with 

M[x2,0}= [ P + 2dx 2 ) =1+ [ ' G(*i,0<k 
Jo Jo 

+ I 2 Q (h, d? P Q (*!,?!) d ft + (30) 
Jo Jo 

where I is the identity operator, and the operator 



M [x 2 , 0] is formally a matricant of ( 28 ) defined in a stan- 
dard fashion through a multiplicative integral J expand- 
ing in the Peano series^. In the same spirit, the operator 
M. [1,0] given by p0| with x% = 1, i.e. taken over a pe- 
riod 1 in x%, may be called a monodromy matrix. It has 
the important property that if e l M". fc i) with k 2 e K is 
an eigenvalue of M. [1,0], then w and k = (ki,k 2 ) satisfy 
Eq. ([5]), i.e. a; 2 is an eigenvalue of ^ with the Flo- 
quet quasi-periodic conditions along both coordinates X\ 
and x 2 . This is similar to the case of scalar waves in 2D 
media with ID periodicity (see^); however, the pres ence 
of terms of the order O {kf) , 0{k\) 3 A in 129} un- 



derlies an essential difference in the 2D periodicity case. 
Note that Ai [a, b] at lj = 0, ki = has the eigenvalue 
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e ik 2 (o,o) _ i corresponding to the eigenvector (1, 0) T , i.e., 
to w (xi) — const. 

The MM approach enables deriving a new form of 
the exact solution for the effective speed. Referring for 
brevity to the isotropic case, it is as follows: 



1 



v (o,i)(^i[i,o]-2:r 1 (i,o) T ;^ 



(31) 



where .Mi [1,0] is .M[1,0] with uj, k x = and 
(.Mi [1, 0] — T) 1 (1, 0) T is any vector from the preimage 
of the vector (1,0) T with respect to Mi [1,0] - 1. We 
will not, however, discuss Eq. (31) in detail because, as 



any exact solution for c, it defies a closed form and hence 
exceeds the scope of the present study. 

Seeking specifically a closed-form estimate of c necessi- 
tates some additional simplifications. On this ground, let 
us further consider the matrix operator A4 which con- 
sists of the first two terms of the Peano series of M. [1, 0] 
(see (30) with x 2 = 1): 



7W[1,0] =M 



with Mo (w, fci) =X+(Q) 



(32) 

Denote by e lk2 and e the eigenvalue and eigenvector of 
.Mo which at uj = 0, k x — coincide with those of 
M [1,0] , so that 



M e (uj, knxi) = e ik ^"' k ^e (uj, k x ; x x ) , 
where k 2 (0, 0) = 0, e (0, 0; x x ) = (1, 0) T . 



(33) 



The motivation for introducing M Q is that k 2 (uj, ki) has 
an exact closed-form asymptotic form that can be used 
for constructing an estimate of c. It is emphasized that 
the difference between M. [1,0] and .Mo, which is given 
by the members of the Peano series ( 30 1 of the order 
n > 2, contains the terms of the same order O (k x ) , 
O (k x ) (9 A) and O (uj 2 ) as in Mq but with numerical 
factors decreasing somewhat like 1/nl. For the latter 
reason, the asymptotics of k 2 (uj, k x ) and k 2 (uj, ki) should 
be close. 

To obtain the asymptotics of k 2 (uj, k x ) in small uj, k x , 
it is convenient to pass from the matrix form of ( 30 1 to 
the scalar equation as follows: 



(Q) X2 e = Ae ( A = e^ - 1 



(34) 



where A = and e\ (x x ) = 1 at uj = 0, k x = by 
(33). Denote fci = ae and uj — fie where e is a small 



perturbation parameter. Inserting in (34 ) 3 and invoking 
p9| yields 



(K + eK x + e 2 TZ 2 ) e x (e, x x ) = A 2 (e) Ve x (e, x x ) 
with Vw = i^ 1 )^ w, TZ w = - ((fi) X2 w')' , 
Tl x w = -ia ({fj) X2 w' + ({fj,) X2 w)'j , 



2 w= ya (v) X2 -/3 (p) X2 )w. 



(35) 



Note that the operators IZi and T> acting on w (x x ) G W 
are self-adjoint with respect to the inner product (/, h) = 
J Q fh*dx x = (fh*) x where * means complex conjugate. 
Applying the standard technique of perturbation theory 
then leads to 

A 2 (e) = (A 2 )i£ + (A 2 ) 2 e 2 + 0(e 3 ) with 
(7?.ie i, e i) 



(A 2 )i = 



(T>e i,e i) 

(TZ 2 e 0X , e i) - (^^leoi^ieoi) 
(£>e i,e i) 



(36) 



where eoi = ei(0, x x ) = 1. First note that (7?-ieoi, eoi) = 
— ia J Q ({p) X2 )' dx x — since ((Jt) x is a periodic function 
of x x , and hence (A 2 ) x = 0. To find TZQ 1 TZ x e 0X = (f (x x ), 
we need to solve the equation lZ (f> (x x ) = lZ x e 0X , that is, 

-^) x J)' = - ia (^)J =, 



(j)(x x )=K + iax x +K x (n) X2 dx x , (37) 
Jo 



where K and K x are constants. Using the boundary con- 
dition 0(0) = (f)(1) for (f>(x x ) € W (see ([27])) determines 
K x , whence 

7?. ~ 1 7?.ieoi = (f ( x i) — K + iotx x 

-ia £ (»}-}dx x . (38) 

Thus, calculating 

(IZq 1Z x eoi,7£ieoi) 
(7?. 2 e i,e i) 
(2?e i, e i) 



(3 9 ) 



and inserting in ( 36 ) yields the explicit form of A 2 (e) = 
( e ik 2 py (A 2 ) 2 e 2 which, on reverting to the original 

parameters k x = ae and uj = fie, yields 



k 2 x "-uj 2 ( P ) 



0(k\,uj^ 



(40) 



As argued above, k 2 (a;, k x ) at small uj and k x is supposed 
to be close to k 2 (uj, k x ); therefore replacing k 2 in (40) 
by k 2 leads to the approximation for the effective speed 
c (k) — lim LJ .fe_j.o w (k) /k (k = fcre) as 

o^ w ^(,;(w-'):; +.!(<,->;:) j. («> 

Note that applying the same scheme with respect to 
the reverse order of coordinates, i.e. imposing the Flo- 
quet condition along x 2 and using the monodromy matrix 



G 



along jci, yields the formula which follows from (41 1 by 
interchanging x\ xi and K\ <^ K2- Neither of the two 
approximations is generally preferable, so it is natural to 



use their average, say, the half-sum | [(41) + (41 (1^2] 



Thereby we arrive at the estimate for the effective speed 
in the following form: 



c 2 (k) 



1 

W) 



(42) 



where (•) is defined by (26 1 (obviously the assumption 
of unit and equal periods Ja^ | is no longer needed) . 
The wave speed estimate (142]) describes an ellipse of ef- 



fective slowness s(k) 



[k)k with the principal axes 



along the translations a! _L a 2 of an orthotropic lat- 
tice. For an isotropic lattice, where each of \{^) x 1 ) 7^ 

(/i -1 ) ) is invariant to x\ ^ X2, Eq. (42) (in con- 



trast to (41)) becomes isotropic, i.e., yields" 
value 



C MM 



1 



2(P> 



(43) 



for any k. It is instructive to apply the explicit formula 
(43) to a square lattice composed of J = 1,2, ... homo- 



geneous materials, which is the case exemplified in detail 
in fjv] Inserting [i (x) = J2j I^jXJ ( x ) f° r x e T, where 
Xj (x) ((xj) = /./) is an indicator function equal to 1 
on the domain occupied by the J th material and to 
elsewhere, reduces Eq. (43) to 



1 



2(P> 



d^ 2 



Ej/Vxjfe) 
_ i 

d<a 

Ej v-jxj (a) 



(44) 



with xjfe) = JoXjfebftOdcL and q = Xi/\&i\. Now 
suppose that one of the constituent materials has jij — > 
and it is distributed with a small (but finite) concentra- 
tion fj along th e un it-cell boundary. Then both integrals 
on the r.h.s. of (44) tend to zero, and so c 2 AM —> 0. Thus 



the essential attribute of the MM estimate (44) is that 



it is capable of capturing the 'insulating' effect of even 
a small concentration of soft material when this forms a 
'network' breaking the connectivity of stiff components 
in the lattice. One more revealing example is the limiting 
case where /i(x) is constant along some fixed direction in 
M 2 (while p(x) may remain 2D-periodic). Taking this di- 
rection as the base vector ei implies (pt) Xl — ^ (2:2) and 
thus reduces ( 42 ) to the well-known exact formula 



(45) 



In fact, the original non-symmetric estimate (41 ) reduces 



to the exact form (45) when /Lt(x) is constant along the 
direction e,-, i = 1 or 2. 



In conclusion, a few remarks are in order concerning 
the approximate nature of the MM-approach implemen- 
tation and result. First, the assumption that /x(x) and 
p(x) are smooth can actually be relaxed to include piece- 
wise continuous functions and hence to apply the ap- 
proximation (42) to composites with inclusions, see fvl 



This is similar to the effect of truncating PWE series of 
piecewise continuous /i(x) and /c(x), which allows one to 
think of them as smooth functions QIV). A second re- 



mark is that the MM estimate (42) is not restricted to 



the isotropic case like the PWE estimate (24>) is. On the 
other hand, due to the simplification adopted on deriving 
Eq. (42 ), it does not contain a cross term proportional 



to K1K2 and hence is unable to pinpoint the effect of 
asymmetric form and/or distribution of inclusions in a 
rectangular lattice that could tilt the principal axes of 
the exact effective-speed curve away from the translation 
vectors ai, a2. For the same reason, Eq. (42) may not 



be invariant with respect to different choices of a unit cell 
in a given lattice. Such deficiency could be rectified by 



the same taking into account the terms of order O (w 2 ) from the 
next (n > 2) Peano-series members, which are discarded 
in Mo (cf. (30) and (32)); however, this is hardly an 



expedient course of action since adding even one more 
term on top of Aio leads to quite a cumbersome expres- 
sion for c. Finally, we note that instead of taking the 
arithmetic mean leading to ( 42 ) , one could have invoked 



another average, e.g., the geometric mean. Its direct use 



1/2 

as [(41 ) x (41 )i^± 2 ] is unreasonable since the resulting 



estimate of squared speed c (k) would no longer be a 
quadratic form in k; however, the geometric mean could 
be applied separately to the coefficients of Kj thus yield- 
ing 



c 2 («) 



1 



-1\ V2 



= c 2 — 

MM' 



(46) 



Comparison of the two MM estimates and Cmm is 



considered in 9V A 1 



IV. PWE NUMERICAL IMPLEMENTATION 

PWE numerical implementation rests on calculation 
of the quantity M{n) = (B _1 d, d) , Eq. (11), which in- 
volves the inverse of the formally infinite matrix B trun- 
cated in the 2D calculations to a finite A^ 2 x A 2 size (A is 
the number of Fourier terms in one coordinate) . Its inver- 
sion takes O (A 8 ) steps. Calculating B _1 d, i.e. solving a 
linear system Bh = d for unknown h by Gauss or similar 
methods, takes O (A 6 ) steps (and needs O (A 4 ) memory 
cells for storing intermediate results). This may also be 
onerous for large enough A. Note also that the case of 
high-contrast lattices with very soft or void components 



7 



needs special care (see e. gP). The difficulty arises due 
to the fact that B is not invertible if /i (fi) = for some 
domain fi of x within the unit cell T. This does not pre- 
clude numerical inversion of truncated B (since a finite- 
size B can no longer possess eigenfunctions with a sup- 
port in fi S T); however, both inversion of B and solving 
Bh = d with zero or small fi (fi) may become tricky be- 
cause taking more elements of B implies a greater impact 
of its small eigenvalues and thus may actually deteriorate 
numerical accuracy. 

In this light, we advocate the method of direct com- 
putation of M (k) via the series expansion (20j) with 
/i = \ (Mmax + Mmin) = M- On fixing the meaning of 
truncated quantities as defined on a iV 2 -dimension sub- 
space l% 2 C I 2 (r\ {0}) spanned by N 2 = (2j + l) 2 vec- 
tors e g = (^gg')g'^o with < |<7j| < 2nj (i = 1,2), the 
explicit expression for computing M (k) is 



MM^/j-iV™ J(-C N2xN 2) n f N 2,f N2 ), (47) 

^ — 4 n— 



N 2 xN 2 



(m) 



C and f, 



N 2 



f have components 



) and (f , e g ) in l 2 N2 . 'Termwise' (by way of stor- 



where C 

ing C"f and calling" on it for C" +1 f = C (C"f)) calcula- 
tion of (47) takes O (miV 4 ) steps, which is notably less 
than O (N°) when N ^> m, 1. The validity of approxi- 
mation ( |47| can be justified on the basis of the sufficient 
condition [T( 



< 1 for convergence of (20 



and on the 

fact that C is close to the operator of multiplication by 
(m( x ) — mo)/mo whence [|C(/Z)[| ~ |/i( x)//I — 1| (see the 
discussion of Eqs. Q, 0\\ in 01. Thus \\C(jl)\\ is 



expected to be less than 1, being close to 1 in the spe- 
cial case where /i is very small in some € T. In the 
former case, fast convergence of (20>) is facilitated by 
the diagonal predominant structure of I + C (see Ap- 
pendix). In the latter case (small ^t(fi)), the fact that 
|f| decreases as g grows large may come into play. How- 
ever, in contrast to the MM approach (see SIIIB), the 



PWE considerations seem unable to explain the very dif- 
ferent effect of this small fi when it occurs either strictly 
inside the unit cell (soft inclusion) or along its boundaries 
(soft matrix). Numerical examples provided in ^Jvjshow 
that Eq. ( 47 ) is not sensitive to fi of an inclusion tending 



to zero, and hence it can be directly applied to comput- 
ing the effective sh ear sp eed in solid/air and solid/fluid 
composites (see e.g! 18 * 23 |), where the solid phase remains 
connected^. The alternative case, in which a very soft 
matrix material forms an 'insulating network', is known 
to be particularly subtle for any PWE-based numerical 
scheme. No wonder that application of Eq. (47) to this 
case requires more numerical effort as detailed in JV] 



Note that taking (20 2) with fi = ~p, which leads to 
the same form (47) for any /z(x), does not at all guar- 
antee the fastest convergence for all fi(x). This is eluci- 
dated in Appendix which contains an example of strict 
and quantitative convergence analysis of the series (20>) 
for a particular family of periodic /i(x). 



V. DISCUSSION AND EXAMPLES 
A. Two-phase lattices 

1. Estimates 

Consider a 2D square lattice which is isotropically com- 
posed of two homogeneous materials J — 1, 2 with con- 
stant pj, \xj and with filling fractions fj (/i +/2 = 1). It 
will also prove useful to introduce the conjugate lattice 
by the following definition: two conjugated binary lat- 
tices are related to one another through the replacement 
J = 1,2 *^ 2,1 (that is, fix, j\ <^ fi%, /a) interchang- 
ing the materials along with their filling fractions. The 
conjugated lattices are referred to below as 1/2 and 2/1 
lattices, with the matrix material put first. Note that the 
exact effective speeds in conjugated lattices are in general 
certainly different, c^/ 2 ) 7^ c (2/i)j except for particular 



symmetric lattice configurations, see iVA2 

The PWE estimate ([24^) of the effective speed c re- 
duces to the form 



4we = y^y (mi/i + M2/2 



/1/2 (Ml - M2) 

Ml + M2 



. (48) 



which, by definition, depends only on the filling fractions 
fj and is not sensitive to the inclusion shape. It is also 



evident that ( 48 ) is invariant under the interchange J = 
1, 2 <=^ 2, 1, i.e., cpwe is the same for the two conjugated 
binary lattices. 

The MM estimate c|j M for the two-phase square lattice 
is given by (44) with J = 1,2. It is not invariant to 
interchanging J = 1, 2 ^ 2, 1, i.e. the effective speed for 
each of the conjugated lattices has its own MM estimate 

C(i/2) ~ 4iM and cp/i) ~ (where c^jj = 

for the symmetric configurations). 

The estimate obtained by means of the multiple- 
scattering theory (MSTP^is, for the 1/2 lattice, 



"(1/2) 



" (1/2)" 


2 = Ml 


( Ml - 


1- A*2 - 


(Mi 


C MST 


(p) 




V[l2 + 


(Mi 



M2)/2 



J = 1 is matrix, J — 2 is inclusion. 



(49) 



Interchanging the indices J = 1,2^2,1 in (49 1 provides 



the MST estimate for the conjugated 2/1 lattice as 

'2/ii - (fix - fj, 2 ) /2 

J (2/l) 
J 



' (2/1)" 


2 


M2 / 


C MST 




(p) \ 



= 2 is matrix, 



2fj, 2 + (mi - Ma) /2 
J = 1 is inclusion. (50) 



The MST estimate defines distinct values of cmst for 
the two conjugated lattices. The choice as to which of 



the MST formulas (49), (50) to apply to a given binary 



lattice depends crucially on the designation of the two 
constituent materials as 'matrix' and 'inclusion'. There 
is no ambiguity for simple configurations where one of the 
materials ('matrix') circumvents the unit-cell boundary 
and the other is enclosed within ('inclusion'). However, 



in the case of a symmetric lattice configuration, for which 



two conjugated lattices are equivalent, Eqs. (49 1 and ( 50 ) 



provide two starkly different MST approximations of a 
single exact value Cn/2) = c (2/i): see further !VA2 

Note that the explicit expressions (|49|, J50|) actually 
have a long record in micromechanics, see^"^. In particu- 
lar, they are the Hashin-Shtrikman bounds (respectively, 
upper and lower at fii > /i 2 or vice versa at Hi < 1^2) ob- 
tained by the variational approach for a binary composite 
of a matrix material J = 1 or 2 with statistically homo- 
geneous inclusions of material J = 2 or 1, see^. The 
relation of these bounds to periodic structures may not 
be generally obvious. At the same time, for the two-phase 
lattices, it is easy to verify explicitly that the PWE esti- 
mate ( 48 ) is always enclosed between ( 49 ) and ( 50 ) , and 



that the upper Hashin-Shtrikman bound is never greater 
than the PWE bound (24[) for the two-phase case; how- 



ever, the same is not always true for the MM estimate 
(44 ) with J = 1,2. One more general result from the the- 
ory of 2D two-phase c ompos ites is noteworthy, which is 
Keller's duality relatiorPSIlZl for the effective shear coeffi- 
cients /z e ff of two reciprocal lattices {(11,(12) and (/i2,/ii) 
obtained from one another by interchanging pi <^ (i 2 
while keeping the concentrations /1 2 intact (cf. the def- 
inition of conjugated lattices). For the isotropic case in 
hand, this relation yields the identity 



(51) 



Among the above-mentioned estimates of c, the MST for- 
mulas ((49]), Q satisfy while the PWE and MM 
approximations J48"t and (44) do not. Note that the 
MM estimate (|46j does satisfy (51); however, the nu- 
merical tests (omitted from the graphical data below to 
avoid its overloading) show that fitting of c by the MM 
estimate cmm given by (44) is always better than by 



'MM 



< Cmm) given by (46 1. The degree to which it is 



better is often quantitative small, but then the departure 
of Cmm from the duality identity (51 ) is equally small. A 
greater accuracy of ( 44 ) than of ( 46 ) extends to the case 
of J > 2, where (46) has no methodological advantage 



of satisfying (51) since the latter is no longer relevant. 



Thus, all in all the MM estimate in the form (44) ap- 



pears to be preferable to (46 1 



2. Examples 

In this subsection, a comparison between the numerical 
evaluation of the effective speed c and its different esti- 
mates is demonstrated for several examples of a square 
lattice of parallel square rods embedded in a matrix and 
oriented at an angle of 0° or 45° to the translation vec- 
tors. Such configurations of phononic crystals have been 
studied, e.g., irP^l. It is clear that the MST estimate 
pHH3 though derived for cylindrical inclusions, should 
be equally viable for square ones since it describes the 
quasistatic limit. If the contrast of matrix and inclusion 



cpb 



c num 

Cmst 
cpwe 
cmm 



CAl 




1 / 



FIG. 1: Effective speed c versus concentration /ai for con- 
jugated Al/Pb and Pb/Al lattices of 0°-oriented rods. The 
numerical curves for both lattices (computed via (47 1 with 
TV = 25 and m = 10), the PWE estimate PH), the MM esti- 



mate ( 44 1 and the MST approximations ( 49 1, ( 50 \ all merge 



at the scale of the plot. 



shear coefficients is relatively low, then so is the differ- 
ence between the two values of the effective speed c for 
the two conjugated lattices. In this case, the PWE, MM 
and MST estimates dish, (1441 and (|49T)-((50l all yield close 



values that provide a good approximation of c in either 
of the conjugated configurations. This is exemplified in 
Fig. 1 for Al and Pb phases with the material constants 
PAl = 2.7, p Pb = 11.6 g/cm 3 and = 26, ^ Pb = 14.9 
GPa. Note that the series ( |47| ) needs only about j ^ 7 
modes (N ~ 15) and m ~ 5 terms for accurate calcula- 
tion of the numerical curve c (/ai) (the larger values of 
N and m indicated in the caption were taken for better 
precision). 

Addressing the high-contrast case, consider two exam- 
ples of binary materials with a 'medium' and 'drastic' 
contrast: one consisting of steel (= St) and epoxy (= Ep), 
and the other of steel and rubber (= R). Their material 
constants are pst — 7.8, ps P — 1-14, prs, = 1.14 g/cm 3 and 
USt = 80, (L Ep = 1-48, (i R = 4- 10~ 5 GPa. The results for 
the St/Ep and Ep/St conjugated lattices of 0°-orientcd 
rods are shown in Fig. 2a, and the results for the St/R 
and R/St lattices are shown in Fig. 3a. It is seen that 
the two numerical curves c (/) , plotted for each conju- 
gated pair as a function of concentration of the same (say, 
softer) material, have quite different trajectories between 
the fixed end points. The physical reason is obvious: the 
effective speed c is indeed strongly affected by a small 
concentration of a highly contrasting component when 
this forms a 'network' breaking up connectivity of the 
volume-dominating component. On the numerical side, 
given the 'medium-contrast' case of steel-epoxy compos- 
ite, Eq. ( 47 ) provides a reasonable approximation of c (/) 
when taken with j — 7 modes (N ~ 15) and m ~ 50 
terms (compare with the above Al-Pb case). About this 
number of modes and terms in Eq. ( 47 1 is also sufficient 
to capture the shape of the curve cjj) for the 'drastic- 
contrast' steel-rubber structure but only if rubber is an 
inclusion located inside the cell. Markedly more numeri- 
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c num 

c MS T (47), (48) 

CPWE 

CMM 

□ e p Est 




c num 

cmst (47), (48) 

CPWE 

CMM 

□ e p Est 



p=0 



FIG. 2: Effective speed (a) for the conjugated St/Ep and 
Ep/St lattices of 0°-oriented rods and (b) for the symmetric 
St/Ep lattice of 45°-rotated rods. Numerical curves c(/e p ) 
are computed via p7| with T V = 841 and m = 150; the PWE 

the MM estimate cmm is 



estimate cpwe is given by p8| 
given by (441; the MST app 
are given by ( 49 1 and ( 50 1 with J 



given by (44 1; the MST approximations c[^ r Ep ' and ci Ep 



St, J 



-MST 

Ep. 



cal effort is required when rubber is the matrix material 
distributed along the unit-cell boundaries - in this case 
no less than j = 12 modes (N ~ 25) and m ~ 150 terms 
in Eq. (47 1 are needed to obtain good accuracy (see § 
4). Note that formally reducing /ir to zero causes no 
discernible changes at the scale of Figs. 3, 4. 

Let us now examine the PWE, MM and MST esti- 
mates of c for the above examples. It is evident that a 
single curve of the PWE estimate, which depends only on 
volume fraction and disregards geometrical details (see 
j V A 1| , cannot fit two markedly different curves of conju- 
gated lattices. As noted in §III| it must be more accurate 
when the stiff component is volumetrically dominant over 
the soft one rather than when the situation is reversed. 
This is what is observed in Figs. 2a and 3a. It is also seen 



c num 

cmst (47), (48) 

CpWE 

cmm 

□ r n st 




c num 

cmst (47), (48) 

CPWE 

cmm 

□ R □ St 



FIG. 3: The same as in Fig. [5] (a,b) but for St/R and R/St 
lattices. Numerical curves c (/r) are computed via ( |47| l with 
TV = 29 and m = 150. Note that cr is not distinguishable 
from at the scale of the vertical axis. 



that the MM and MST estimates provide a fairly close 
evaluation of c, which fits very well the whole numeri- 
cal curve of c for St/Ep and St/R lattices (soft rods in 
stiff matrix); however, they lose accuracy for the conju- 
gated, Ep/St and R/St lattices (stiff rods in soft matrix), 
specifically when the rod concentration /st (= 1 — /ep.r) 
is close to 1. Regarding MST, this is in agreement with 
the remark made on its derivation irPHED that the MST 
estimate does not fully account for the multiple interac- 
tions and hence may be error prone in the case of densely 
packed stiff inclusions. Thus, in the latter case, the PWE 
estimate is preferable to two others, as illustrated in Fig. 
2a and especially in Fig. 3a. 

Consider next similar structures but with 45°-rotated 
rods, which is the case where the two conjugated lat- 
tices coincide into one symmetric configuration. The 
corresponding dependence of the effective speed versus 
concentration c (/) has a single- valued approximation for 
each of the PWE and MM estimates, whereas the MST 
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cmm c mst 
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FIG. 4: Effective speed as a function of concentration of inclu- 
sions in (a) St/R and (b) R/St conjugated lattices of circular 
cylinders in a mat rix. Numerical curves c(/r) and c(/st) are 
computed via (47 1 with N — 29 and m = 150. 



B. Three-phase lattices 

1. Estimates 

Consider a 2D square lattice similar to above but with 
a coated inclusion. Such nested structures have received 
much attention lately in rel ation to modelling locally res- 
onant phononic crystals, e.g! 24 | 25 l The PWE and MM es- 
timates of the effective speed c for this case are given by 
Eqs. ([24}) and ([44]) with (•) = J2j 0)j fj and J = h 2, 3. 



If the concentration fj of one of the constituent mate- 
rials tends to zero, the MM estimate (44) for three con- 



stituents certainly tends to that for two remaining con 
stituents; whereas the PWE estimate (24 
fs tends to its form for the pair J 



with, 



say, 



1,2 only if the 
'vanishing' material is neither the stiff est nor the softest 
one, i.e. if ^ 3 ^ 

Mmin ; Mmax • 

As a MST counterpart, we adopt the generalization of 



( 49 1 that is well-known in micromechanics as the Kuster- 
Toksoz formula (closely related to Hashin-Shtrikman 
bounds) for 2D flu ids with small concentration of differ- 
ent inclusion d 26 ! 27 !. More recently, it was used for a peri- 
odic structure of different cylinders in a fluid matrisP^] 
The formula for the 2D configurations considered here is 



L— M.J 



2 _ Ml ( 1 J2,J = 2 fj^ 1 1+l : : \ 

1 is matrix, J — 2, 3 are inclusions. 



(52) 



The MST estimate ( 52 1 coincides with the binary formula 



( 49 1 if any one of the inclusion concentrations f 2 or / 3 
is zero. On the other hand, (52) does not tend to either 



of ( 49 ) and ( 50 ) as the matrix concentration f\ tends to 



zero (which is not surprising since the Kuster-Toksoz is 
not recommended at low matrix concentratiorpSI) . 



estimate still defines two different approximations (49) 



and ( 50 1 for the single curve c (/) . Comparing these es- 



timates displayed alongside the numerical curve c (/) in 
Figs. 2b and 3b shows that the PWE estimate is the 
most accurate so long as the stiff component is volume- 
dominant; the MM estimate provides the best 'overall' 
fit; and each of the MST approximations works over less 
than a half of the range while mismatching markedly the 
other half. 



Finally, we consider the case of cylindrical inclusions. 
Results for the steel - rubber conjugate lattices with cir- 
cular rods are presented in Fig. 4. It is instructive to 
observe the similarity of the dependences c (/) on the 
concentration of inclusions / = fst and /r, which are 
displayed in Figs. 4a and 4b, to the two corresponding 
'halves' of the corresponding curves for square rods in 
Fig. 3b. 



2. Examples 

Denote the filling fraction of a coated inclusion in a 
matrix (J = 1) by / and set the filling fractions of the 
skin (J = 2) and core (J = 3) materials as 

f 2 = af (skin), / 3 = (1 - a) f (core) =>> f 2 + f 3 = f. 

(53) 

The effective speed c of the three-phase composite is now 
a function of the single variable / = 1 — f\. 

Motivated bjEES we first examine the case of a soft 
coating (skin) material. Consider the square St/R/Pb 
lattice of square lead (=Pb) rods coated by rubber (=R) 
which are embedded in steel matrix (Fig. 5a) . The value 
of c (/) at / = is obviously the speed in the matrix, 
c(0) = est- The opposite limit value of c(/) at / = 1 
is equal to the effective speed CR/pb (/r) in the binary 
R/Pb lattice of lead rods embedded in the rubber matrix 
with the volume fractions fixed by ( |53[ ) as /p, = a and 
/pb = 1 — a. Once /r is not too small, c R / Pb (/r) should 
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FIG. 5: Effective speed c (/) for three-phase lattices where 
/ is given by Q with a = 4/9: (a) Pb/Ru/St structure of 
coated square rods and (b) Ep/St/Ep structure of cylindrical 
annuli. Numerical curves are computed via ( |47[ ) with TV = 29 
and m = 150 



be close to cr (see Fig. 3a) , which therefore implies that 
c(/) in the St/R/Pb structure has a very small value in 
the limit / — s- 1. This is observed in Fig. 5a (where a = 
4/9). It is also seen that the PWE and MST estimates 



(24>) and (52) of c(/) do not describe this behaviour of 
c[7) at / — > 1 and overestimate c(l) (by an incidentally 
close value which is neither PWE nor MST estimate of 

above) 



CR/Pb (/r), as pointed out in S 



the MM estimate ( 44 1 provides a g< 



VB1 



By contrast, 
food fit for the whole 
curve c (/) including the critical region / — > 1. This 
is because Eq. (44 1 captures the 'insulating' effect of 



a small concentration of soft material which drastically 
decreases the effective speed when this material extends 



throughout the unit-cell boundary, see [MB 



Another case of interest is when the matrix material 
coincides with that of the rod core, which means that 
the rod coatings are simply spacers separating the same 
material. Figure 5b demonstrates the dependence of the 
effective speed c on the concentration / of stiff (steel) 
cylindrical annuli embedded in a soft (epoxy=Ep) mate- 
rial. The shape of the curve c (/) can be shown to change 
only slightly if the steel spacers are square instead of cir- 
cular. It is seen from Fig. 5b that the basic outline of this 



curve is again best approximated by the MM estimate. 



VI. CONCLUSION 

The paper uses the PWE approach and a newly devel- 
oped MM approach, based on the monodromy matrix, 
to derive the new estimates of the effective shear-wave 
speed c in 2D periodic lattices. The estimates are com- 
pared with the known MST approximations and with the 
numerical data for a number of examples of two- and 
three-phase square lattices. The main findings are listed 
in the Introduction. The results for effective velocities of 
the vector waves in the 3D lattices are to be reported else- 
where. It is worth pointing out that the obtained PWE 
and MM estimates are also valid for the gradient-index, 
or functionally graded, materials (for which the MST is 
irrelevant). In conclusion, the combination of the pertur- 
bation theory with the PWE and MM techniques, which 
is elaborated in this paper, is hoped to lend an efficient 
tool for a broad range of problems concerned with peri- 
odic composites, phononic crystals and metamaterials. 
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APPENDIX. Convergence of (|20j>): a strict ex- 
ample 

Sufficient condition on /i(x). Our objective is to 
provide a rigorous example of a class of functions p, (x) = 
Mo + (x) that guarantee convergence for M(n) — 
Mo^S^Lo (( — C) n f, f), and thus validate application of 
this series for computing the effective parameters /i e ff (k) 
and c 2 (k). To do so, we begin by formulating a suffi- 
cient condition on (i (x) to fulfill the sufficient condition 
||C(/xo)|| < 1 for convergence of (pO^) as m — > oo. Note 
that the matrix C can be written as 



Mo 



Jgfe.g'] 



Ma (g) Jg where J g : 

if g = g - g', 
otherwise. 




(54) 



It is seen from (54) that ||C|| < fi Q X) g er I MA (g)| since 
||jg|| < 1, which in turn is because all its nonzero ele- 
ments Jg [g, g'] occupy a single particular diagonal and 
satisfy |Jg[g, g'] < 1. Hence the sufficient convergence 
condition ||C|| < 1 may be eased to 



|C( W )||<) Mo" 1 £ ser lMA (81 = 1 



1 E 



<M> 
Mo 



|/i(g)|s6 w <lfor mo>(m>. (55) 



^0 ^8^0 

In other words, for those /x(x) which satisfy 



£ |£(g)|<<M> 



(56) 
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there always exists a choice of > (m) which ensures 
||C(^o)|| < 1 and hence guarantees convergence of (20>) 



to M(k). The remainder of the series (p0|) with /i > 
(ji) may be estimated as follows 



Mo _1 E ,,((-C) B f,f) 

' * 7 ■ 1 



< 



'n— m+1 

<m> 2 ©r 1 



< 



Mo 
<M> 2 © 



E 



n—m+l 



2 *n»m+l 

Mo 



mo i - ©mo (m) - Eg^o Im(s) 



(57) 



where it has been used that ||C|| < Mo by (551 and that 



|/(x)| 2 )<max|/(x) 



(58) 



< 



E 



<<M> 



for /(x) = E g ^o/(g)e ig x by_Jl4( and 
least value of the residual sum (|57| for all /zq 



(56) 



The 

> (n) is 



achieved when Mo is minimum, which is the case when 
Mo = (m)- Note that the average (fj) of /i(x) satisfying 



(561 may well differ (be greater or less) th an t he value 



as a nu- 



M = i (Mmax + Mmin) , which was argued in fIV 
merically reliable choice of (j,q in (pO^). There is indeed 
no contradiction in this difference. First, recall that all 
the conclusions of Appendix stem from only the sufficient 
conditions. Second, as mentioned in |IV| an advantage 
of taking (|20^ ) with /j,q = fL is that it yields the same 
formula ( |47[ ) for any profile m( x )i but this choice of /j,q 
is not intended to provide the fastest convergence for all 
possible profiles. 

We still need to examine the restrictions on /i (x) 
whi ch are imposed by the derived sufficient condition 
d56b. First of all, by g /x(x) = £ g M(g)e igx > 
mTO) - Egio lM(g)l > °i Le - onl y positive /z(x) are al- 



lowed as needed. Second, any ^t(x) satisfying (56) must 
have a uniformly converging Fourier series and hence be 
continuous. The latter is actually not a loss of generality 
in the numerical context, even if we are mostly interested 
in the case of materials with inclusions (i.e. with jumps of 
properties) , because the calculations deal with truncated 



Fourier series of A*(x) which in effect replaces a possibly 
piecewise constant /i(x) by a continuous profile. Thirdly, 
(561 implies that |^t(x) — (fi)\ < (fi), i.e. /i(x) > should 
not depart 'too far' from its average (/i) . When so, the 
matrix I + C is diagonal predominant and |f | decreases 
for large g, both furthering the truncation of the PWE 
and of the power series in (20;). It is evident that the 



above condition, which may be recast as ^t max < 2 (/z) , 
fits a fairly broad class of functions /k(x). 

Example. In constructing an explicit example of the 



profile /i(x) which ensures convergence of (20>), we con 



sider one that emulates a high-contrast composite with 
a small volume fraction of soft inclusions. For brevity 
of writing, let T = [— tt,7t] so that x = (2:1, £2), 
g = (91,92) with Xi e [-7r,7r] and g { = n t (rij € Z). 
Denote 

iPmn 2 ( x ) = {xi)i>n 2 (x 2 ) = 

^ — '|si|<»i; |92|<"2 

where ^» (a:) = (cosx) 2 " = V $ n (g) e l9X . (59) 

L — 'lsl<« 

Since ip n (ttI) = 1 (i £ Z) and ip n (x) — > for x 7^ ttI 
as n — > 00, the function tp nin2 (x\,X2) for large n\, n 2 
tends to a 2D grid of narrow unit peaks. Note also that 
$n (3) > V5 and E| g |< n ^n (3) = V'n (0) = 1, whence 
0n in2 (g) > and E| ffl |< ni , \ g2 \<n 2 (g) = L Using 

this (Pmn 2 , define the function /^(x) as follows: 

/x(x) = fi Q + /i A (x) : 

fi > A> 0, ^ A (x) = -Aip nin2 (x) , (60) 

where /io and ^4 are some constants. From the above 
properties it follows that 



V | M (g)| = AV £ ril „ 2 (g) 



A (1 - ^„ in2 (0)) < Mo ~ ^mn 2 (0) = (/x) 



Thus the function (60) satisfies the condition (56) suffi 



cient for convergence of (20>) 
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Generally the condition for such isotropic behavior at any 
d > 2 may be stated as invariance of the coefficient ^i(x) 
in |5| to the shift x\ — > X2, --,Xd — > x\ and, separately, to 
the change of sign x± — > —xi of the Cartesian coordinates 

i,ofxel J . 

Note that the effective density for shear (SH) waves in 
2D solid-fluid structures depends only on the solid density 
since the vanishing shear force on the fluid/solid interface 
means that the fluid does not participate in the SH motion. 



